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ABSTRACT 

We analyze the angular clustering of z 2.3 distant red galaxies (DRGs) measured by lOuadri et al.l (l2008l) . 
We find that, with robust estimates of the measurement errors and realistic halo occupation distribution mod- 
eling, the measured clustering can be well fit within standard halo occupation models, in contrast to previous 
results. However, in order to fit the strong break in w{9) at 6* = 10", nearly all satellite galaxies in the DRG 
luminosity range are required to be DRGs. Within this luminosity-threshold sample, the fraction of galaxies 
that are DRGs is ^ 44%, implying that the formation of DRGs is more efficient for satellite galaxies than for 
central galaxies. Despite the evolved stellar populations contained within DRGs at z = 2.3, 90% of satellite 
galaxies in the DRG luminosity range have been accreted within 500 Myr. Thus, satellite DRGs must have 
known they would become satellites well before the time of their accretion. This implies that the formation of 
DRGs correlates with large-scale environment at fixed halo mass, although the large-scale bias of DRGs can be 
well fit without such assumptions. Further data are required to resolve this issue. Using the observational esti- 
mate that ~ 30% of DRGs have no ongoing star formation, we infer a timescale for star formation quenching 
for satellite galaxies of 450 Myr, although the uncertainty on this number is large. However, unless all non-star 
forming satellite DRGs were quenched before accretion, the quenching timescale is significantly shorter than 
z ~ estimates. Down to the completeness limit of the Quadri et al sample, we find that the halo masses of 
central DRGs are ~ 50% higher than non-DRGs in the same luminosity range, but at the highest halo masses 
the central galaxies are DRGs only ~ 2/3 of the time. 

Subject headings: cosmology: theory — galaxiesxlustering — galaxies: halos — galaxies:formation — large- 
scale structure of universe 



^ , 1. INTRODUCTION 

OO ' Distant red galaxies (hereafter DRGs) are bright near- 
t:J" , inf rared-selected galaxies that have been detect e d out to z ^ 
■ 3 (iFranx et alj 120031: Ivan Dokkum etal] 120031: lOaddi et alj 
l2003r These galaxies are massive o bjects containing a sig- 
, nificant population of evolved stars (iForster Schreiber et alJ 
O ■ I2OOI . with many galaxies hav ing highly suppressed star for- 
! mation or possibly none at all jLabbe et alj|2005l: iKriek et al] 
. I2OO6I) . Given the remarkable fact that these evolved (and pos- 
sibly dead) galaxies existed a mere 2 Gyr after the big bang, 
many studies have further sought to measure their abundance 
I and clust ering to determine their connec tion to dark matter 
' structure (iDaddi et al.ll2003HGr azian et ani2006l:IOuadri et alJ 
C3 ■ [2007', '2008': 'Marche sini etalEoOTtilchikawa et alj|20'07ll 

Daddi et al. (2003^ measured a very strong clustering 
strength for DRGs, implying that these galaxies live in only 
the most massive dark matter halos that exist at such high red- 
shift. In fact, the number density of such halos was signifi- 
cantly smal ler than the observe d space density of the DRGs 
themselves. iDaddi et"an (l2003h used measurements of small- 
scale clustering to find the correlation length of their sam- 
ple through extrapolation to larger scales. Using the halo 
occupation distribution to model the cluster i ng (H OD; see 
ICoorav&Shethll2002l for a review), IZhend ( |2004|) demon- 
strated that the high clustering strength and observed number 
density were compatible by taking into account a small frac- 
tion of DRGs that are satellite galaxies orbiting around other 
DRGs that sit at the center of the dark matter halo. This small 
satellite fraction produces a significant difference in the slope 
of the correlation function at large and small scales. 



lOuadri et all (l2008l) (hereafter Q08) used a much larger 
sample of DRGs to confirm the strong clustering of 
iDaddi et"an (l2003l) and subsequent studies, both at large (~ 20 
/z~^Mpc) and small (^ 200 /!"'kpc) scales. They attempted 
to m odel the data w ith a somewhat different HOD approach 
from iZhengj (|2004|) . but they could not find a satisfactory 
fit: models that fit at small scales underpredicted the large- 
scale bias, while models that reproduced the large-scale bias 
overpredicted the small-scale clustering. They proposed that 
standard assumptions of HOD models do not describe these 
galaxies; their halo occupation correlates with large-scale en- 
viron ment as well as halo mass, i.e., the so-called "assembly 
bias" (ICroton et al.ll2007l:rWechsler et alJl2006l: iGao & Whitd 
l2007h . This correlation would imply that galaxy formation 
depends on large-scale environment a t fixed halo rnass, a cor- 
relation not seen in low-redshift data ( iTinker et alj|2008bl) . 

In this paper we demonstrate that the observations of Q08 
can be well-fit by standard HOD models. The inclusion of 
cosmic variance in the large-scale clustering errors amelio- 
rates much of the discrepancy with their model fits. We also 
show that a more physically motivated halo occupation model 
yields significantly better results in fitting the data as well. 
Further, we will use our HOD modeling to test the mecha- 
nisms for the formation of red galaxies at z > 2. 

Unless otherwise stated, all calculations adopt a flat ACDM 
cosmology c onsistent with the la test constraints from CMB 
anisotropics (Du nklev et alj|2008l) . Our cosmological param- 
eter set is (a,„crg,/2,«„f7i) = (0.25,0.8,0.7,0.95,0.045). AH 
distances are comoving. Due to their color selection, through- 
out this paper we will use the terms DRGs and "red galaxies" 
interchangeably. 
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Fig. 1. — Panel (a): Halo occupation functions produced by the abundance matching method described in the text. Gray symbols represent all halos 
(parent+sub) while white symbols represent only the subhalos. The two simulations have two different cosmologies, so they have been shifted to a common mass 
scale to demonstrate self- similarity. The solid and dotted curves represent the HOD used in all analytic calculations henceforth. Panel (b): Comparison between 
the clustering in the L120 box and the analytic model. The HOD for the analytic model is obtained from fitting the simulation results (e.g., panel a, but shifted to 
the proper mass scale). The solid line is the full coiTelation function, and the dotted hne the one-halo term. 



2. METHODS 

2.1. Definition of the Term "Halo" 

Because we will be dealing with halos that exist in various 
environments, it is important to have a clear definition of a 
halo and discuss what it implies. We assume that all galaxies 
live at the center of a virialized clump of dark matter That 
dark matter clump may be isolated or it may exist within the 
virial radius of a larger structure. Therefore we will use the 
term halo to refer to an object that is distinct; i.e., it does not 
exist within the virial radius of another object. These objects 
typically have a mean overdensity of ^200 times that of the 
background universe. We refer to objects inside the virial ra- 
dius of halos as subhalos. We use the term galactic halos to 
refer to all halos, both halos and subhalos, that contain galax- 
ies at their center. 

2.2. Halo Occupation from Simulations 

Although collisionless A^-body simulations do not include 
any baryon physics, one can associate the likely sites of 
galaxy formation with the dark matter halos and subhalos 
within a simulation. Several recent studies have demonstrated 
the robustness of this assumption by comparing the cluster- 
ing of galaxies to that of a sample of galactic halos with the 
same space density; i.e., galaxies brighter than a given lumi- 
nosity threshold compared with galactic halos more massive 
than a threshold that yields the same abundance. IConrov et alj 
(120061) found that the predicted galactic halo clustering was 
consistent with galaxy tw o-point clustering measurements 
from z = to z = 5 (see als o lKravtsov et al.|[2004l ; FWing et alj 
120061) . iMarm et"an (120081) extended this to measurements of 
the galaxy three-point correlation function as well. 

For the purpose of this paper, we use high-resolution cos- 
mological A^-body simulations to guide our choice of the halo 
occupation of all galaxies (DRGs and non-DRGs). The HOD 
for all galaxies is then fixed by making use of the luminosity 
function of all galaxies (see below), and we focus our effort 



TABLE 1 
List of Simulations 



Lbox (h ' Mpc) 


(r2m,t78,f!j) 


nip [h ' M0] 


Zout 


120 


(0.3,0. 


9,1.0) 


1.07 X 10' 


2.0 


160 


(0.24,0.75,0.95) 


2.54 X 10* 


2.5 


1000 


(0.27,0. 


8,0.95) 


6.98 X 10'" 


2.5 



Note. — Each simulation will be referred to in the text by its box 
size. A ll simulations were performed with the ART code of lKravtsov et al.l 
119971). The L12 Q and LI 000 simulations have been described in 
ITinker et alj )2Q08al) . 



on constraining the DRG HOD from the clustering data. 

Just knowing where the galaxies are, however, doesn't 
identify which ones are red. Before creating a model for 
the halo occupation of DRGs, we first use the subhalo 
abundance matching technique (SHAM) to set the occu- 
pation of all galaxies down to the completeness limit of 
Q08. The space density of DRGs in the Q08 sample is 
«DRG = 6.5 X 10"^ (/?"' Mpc)""* down to their completeness 
li mit of < 21. Using the z ^ 23 luminosity function 
of iMarchesini et al] (|2007|) . the space density of all galaxies 
is 1.5 X 10"^ (/j~' Mpc)"-' at the same magnitude threshold, 
Mr = -22.3, yielding a DRG fraction of 44%.' Whenever re- 
ferring to the sample of all galaxies, we mean all galaxies 
(DRGs and non-DRGs) down to the completeness limit of the 
Q08 sample. 

Figure [T^ shows the halo occupation functions, (A^) m, of 

* A s Q08 point out, the space density of DRGs from IMarchesini et alj 
)2007l) is slightly lower than that of the larger Q08 sample. Thus, whenever 
using the lMarchesini et all )2007ft luminosity functions and data, a correction 
factor of 6.5/5 ~ 1 .3 is applied. This in creases the p ublished number density 
of galaxies brighter than Mr = -22.3 in lMarchesini et ah t2007i) from 1.2 X 
10-3 to 1 5 X 10-3 (/j-l ]y[p(.)-3 
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halos and subhalos above a given maximum circular veloc- 
ity such that the space density of all galaxies down to the Q08 
limit is obtained. We use results from the first two simulations 
listed in Table 1 . Halos and subhalos are identified in the si m- 
ulation by the algorithm described in lKravtsov et al.l (l2004l) . a 
variant of the spherical overdensity halo finder. Because the 
simulations have different cosmologies and redshifts, {N)m 
will differ between the two. However, the shape of (A^) m is 
self-similar. The HOD from both simulations is well approx- 
imated by a satellite occupation function of the form 



(Mat)M=( — 



expl — —l^ (1) 
where asat = 1 , and a central occupation function of the form 



1 if M > M„. 
if M < M„. 



(2) 



Both simulations exhibit the mass ratios Mi /M^in =15.7 and 
Mcut/M,^in = 1.14. To match the number density of all galax- 
ies in our fiducial cosmology, a value of M^in = 9 x 10" 
/i"'M0 is required. This yields values of Mi = 1.4 x 10'^^ 
Mq, and Mcut = 1 .0 x 10'^ /r' Mq. In Figure[IJi, the HOD 
from LI 20, a simulation with WMAPl cosmological param- 
eters dSpergel et al.l l2003), has been shifted by -0.25 dex in 
halo mass. The results from L160, a si mulation closer to th e 
WMAP3 cosmological parameter set (ISpergel et al.l l2007h . 
have been shifted by +0.23 dex. These shifts bring the HODs 
into alignment with the HOD for our fiducial cosmology. In 
the L120 box at z = 2, the total fraction of galactic halos 
that are sub halos is 13%. This is in good agreement with 
IConrov et al . (2006 ) but some what smaller than the subhalo 
fraction of IWetzel et~aLl (l2008h , who find 18%. This is most 
likely attributable to differing classifications of a subhalo. In 
our definition, a halo becomes a subhal o when it passes th e 
spherical virial radius of a larger halo. In lWetzel et al.l (120081) . 
halos become subhalos when they are linked by the friends- 
of-friends algorithm, which can link objects outside of what 
we have defined as the virial radius. 

With {N)m specified, we can analytically calculate the 
galaxy two-point autocorrelation ^(r). O ur mode l for f jr) is 
based on th e analy tic model detailed in IZhend (120041) and 
iTinker et al.l (120051) . We use a Poisson distribution for the 
scatter about (Ns at)M, in agreement with simulations and semi- 
analytic models dKravtsov et al .' 2004; Zheng et al. 2005) and 
with observatio ns of galaxy clusters (Kochaneket al. 2003c 
iLin et al.ll2004l) . We also assume that the satellite galaxies 
follow the radial distribut ion of dark matter w ithin the halos 
(e.g., the density p rofile oflNavarro et al.l 19961) . u sing the halo 
concen trations of IZhao et al.l 120081 We use the iTinker et al.l 
(l2008i) mass function for A = 200 and a new halo bias func- 
tion based on those simulations (Tinker et. al., in prep; see Ap- 
pendix A). Our analytic model incorporates scale-dependent 
halo bias and halo exclusion for proper modeling of the tran- 
sition between the one-halo term and the two-halo ter m. It has 
been fully tested against numerical s imulations dZehavi et al.l 
l2004t ITinker et al.ll200llChenl2007l; Wechsler et. al. 2009, in 
preparation). 

The clustering statistic measured by Q08 is the angular 
correlation function w{d). The three-dimensional correlation 
function is connected to the angular clustering by 



w(e) 



= J dzN\z) ^Jdx£. (Vx2 + r26'2) , (3) 



where N{z) is the normalized redshift distribution of the 
galaxy sample, r is the comoving radial distance at red- 
shift z and dr/dz = {c / Ho) / y^njT+fp + n^. The photo- 
metric redshift distribution is roughly a top-hat function for 
2.0 < z < 2.8. However, photometric redshift uncertainties 
are important in estimating the true underlying redshift dis- 
tribution. Unless specified, we use the estimated true redshift 
distribution of the sample used by Q08 (kindly provided by R. 
Quadri). For this N{z), the redshift distribution peaks at z = 2.3 
but contains significant wings out to z ~ 1.5 and z 3.3 (see 
Figure 3 in Q08 for comparison). The amplitude of the w{6) 
is somewhat sensitive to the choice of A^(z); for the same ^(r), 
using the photometric N{z) increases the amplitude of w{6) 
by ^ 20%, independent of scale. Although the uncertainty 
in N{z) is a source of error in the theoretical calculations, we 
will show that it is subdominant to the sample variance of 
the observations themselves. In practice, we calculate ^(r) 
at the peak of N{z), z = 2.3, and implement that function in 
equation (O. Q08 recalculate ^(r) at each z using the redshift- 
dependent halo statistics with a fixed HOD. The assumption 
that the HOD is fixed in redshift is likely wrong to some de- 
gree but cannot be quantified, thus it is not necessarily more 
robust to recalculate ^(r) as a function of redshift. In practice, 
we find that there are negligible differences between the two 
approaches. 

Our motivation for using an analytic model to calculate 
galaxy clustering, rather than using one of the simulations, 
is two-fold. First, it frees us to use our desired cosmological 
model. Second, it ameliorates any numerical issues of spatial 
resolution; the Q08 clustering measurements extend down to 
^ 0.02 /i"' Mpc (comoving at z = 2.3), which is below the limit 
at which subhalos can be identified in the simulations. Fig- 
ure [T]5 demonstrates the veracity of our analytic model. The 
points show the clustering of subhalos in the L120 box while 
the solid curve shows the clustering obtained from the ana- 
lytic model using the HOD in Figure[T^ (shifted to the proper 
mass scale, and using the cosmology of the L120 box). The 
good agreement with the simulation results demonstrate the 
robustness of our analytic model. 

2.3. Breaking Galaxies into Red and Blue 

DRGs do not comprise a lumino sity-threshold sarnple of 
galaxies; the luminosity functions of iMarchesini et aP (l2007l) 
demonstrate that even at the brightest end, DRGs only account 
for roughly half of all galaxies. Our model for the fraction of 
central galaxies that are DRGs has a mass-dependence of the 
form 



/«cen(M) = /Rn,axeXp 



-kM„ 



M-M„ 



(4) 



Equation (|4|i implies that all galaxies at M = M^m are blue 
(non-DRGs), and that the fraction of DRGs smoothly in- 
creases to an asymptotic value of /smax at high masses. If 
K = 0, then /ficen = /smax at all masses and central DRGs are 
a random subsample of all central galaxies. Formally, the pa- 
rameter K, is the fractional increase in halo mass with respect 
to Mmin for /scen to rcach 1 /e of the asymptotic value /fimax- In 
lieu of K, we will refer to a physically more interesting quan- 
tity, the ratio between the mean halo mass for central DRGs 
and the mean halo mass for all central galaxies, 



Ate, 



X^_^ dM(dn/dM)MfR,,,{M) 
dM(dn/dM)M ' 



(5) 
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Fig. 2. — Dispersion of the clustering from the 15 mock catalogs. The 
thin gray curves show w(8) for a random 18% of the mocks. The open circles 
are the Q08 data and errors. The open squares represent the mean (connected 
for the integral constraint) and dispersion among the mocks. The squares 
have been shifted slightly for clarity. For all calculations, the covariance of 
the mocks is used as the errors on the Q08 data. 



Fig. 3. — Comparison of models to n'(6) measurements. The solid curve 
shows the best-fit HOD model that matches the w(8) and number density of 
DRGs. The dashed and dashed-dot curves show the breakdown of the one- 
halo term into central-satellite pairs (^jj,) and satellite-satellite pairs 
The thin dotted curves ai'e the two fits presented in Q08. 



where dn /dM is the halo mass function. 

For satellite galaxies we stipulate that a constant fraction 
/fisat of satellites are DRGs, independent of mass. A more 
physical model may have a mass-dependent satellite red frac- 
tion, but with the large errors on the given data, a three- 
parameter model is a reasonable first step. 

An important uncertainty in modeling DRG clustering is 
the uncertainty in the number density of such objects. Our 
model for the DRG fraction has three free parameters (k, 
/fimax, and /fisat), but «DRG rcduces the degrees of freedom by 
one. However, given the small sample size, the error in /Idrg 
is non-negligible and must be taken into account when mod- 
eling the data. We will address this in the following section. 

2.4. Error Estimation 

Q08 obtain error bars through bootstrap resampling of their 
data, where each bootstrap sample draws from the distribution 
of all galaxies. This method is good for estimating errors due 
to shot noise but not due to cosmic variance. The measure- 
ments of w{9) are taken out to ^ 1200 arcsec (0.33 deg), close 
to half the width of the field (0.84 deg). Given the volume 
of the field, roughly 2.2 x 10'' (/2"'Mpc)'' (about the volume 
of a (130 /!~'Mpc)-' cube), cosmic variance is important for 
modeling the data at large-scales. 

To estimate the errors on the Q08 data, we construct mock 
catalogs by populating the halos in a large-volume simulation 
(LIOOO listed in Table 1). The simulation itself does not re- 
solve substructure, so we use the analytic (A^) m from Figure 
[T]to populate the halos in the simulation. Central galaxies are 
placed at the center of mass of the halo and satellite galaxies 
are pl aced randomly, assum ing a density profile that follows 
that of iNavarro et al.l (Il996l) . Halos are identified in the simu- 
lation using the friends-of-friends technique (e.g.. lDavis et al.l 
Il985h . Due to the low mass resolution of the simulation, we 
use halos down to 10 particles. This is not ideal, but the 
large-scale bias of these halos matches theoretical expecta- 
tions, thus the mock-to-mock variations of these halos is accu- 



rate jHu & Kravtsovll2003h . To determine the DRG fraction, 
we use a model with parameters k = 0.65 and /^sat = 0.64, 
yielding a value of fumax = 0.95. (This was the best-fit model 
from fitting the Q08 data using the Q08 error bars and a con- 
stant N{z) from 2 <z< 2.8.) 

We use the simulation output at z = 2.5.^ The comoving dis- 
tance from 2 < z < 3 is roughly 1 comoving Gpc, allowing 
us to fully incorporate the the redshift depth of the Q08 photo- 
metric sample. We use a redshift distribution that is constant 
from z = 2 to z = 2.8, and zero at higher and lower redshift. At 
z = 2.8, 0.84deg of arc is 64 comoving /r' Mpc, thus we are 
able to create 15^ = 225 mocks. 

A random sample of 18% of the mocks are compared to 
the Q08 data in Figured The data has been corrected for the 
integral constraint. We also perform the integral constraint 
correction for the correlation functions from the mocks by 
adding a constant. The constant is computed as the differ- 
ence between the mean correlation function of all mocks and 
the correlation function measured from the entire box. The 
agreement between the mocks and the Q08 data at large scales 
is artificially enhanced due to a number of reasons. 

The narrow N{z) function increases the amplitude of w{9), 
and with the smaller error bars at large scales the best-fit 
model has a lower number density, yielding a higher bias. 
The FOF halos with low particle numbers e xhibit some- 
what s tronger scale-dependence than that seen in lTinker et alj 
(I2009I) . which helps the fit at 6* 20". However, for the pur- 
pose of estimating cosmic variance the most important cri- 
terion is producing mocks the reproduce the data on most 
scales. The large-scale clustering of the mock DRG sam- 
ples is consistent with the data given the sample variance at 
less than I -a. At 6* < 50 arcsec, the bootstrap errors of Q08 

^ Using a single redshift output rather than a full lightcone makes neghgi- 
ble difference in the clustering. As redshift decreases, the amplitude of dark 
matter clustering increases with the growth factor, but if the HOD is fixed (as 
assumed here and in all HOD fits to w(6) data) then the bias at the mean mass 
scale of the galaxies decreases in such a way as to nearly cancel the change 
in dark matter clustering. 



are quite accurate in recovering the error, but the error bar on 
the 9 = 10^ arcsec datum is underestimated by a factor of 5. 
Henceforth, we use the dispersion among the mocks as the di- 
agonal errors on the Q08 data. We note that we use the disper- 
sion as an absolute error, not a fractional error. Because w{9) 
is a projected quantity, the data points are highly correlated. 
Our large number of mocks allows us to robustly calculate the 
covariance matrix for the Q08 sample. For all model fitting 
we use the covariance matrix. We also note that Q08 have an 
additional datum ?X 6 = 1400 arcsec. Although the measured 
value of w(9) is significantly lower than the other data points, 
we do not use this point because sample variance renders it 
e ssentially usele ss. 

iGuo & White! (l2008h use the Millennium semi-analytic 
model to make similar mock samples to compare to the Q08 
data (a simulation 1 /8 the volume of our LIOOO simulation). 
The dispersion of their 48 mocks is comparable to Figure |2] 
but their mean correlation function is significantly lower than 
the Q08 data, even at small scales. Even with this discrep- 
ancy, 6 of their 48 mocks match the large-scale data up to the 
data aiO < lO-', which is roughly a 1-cr result and consistent 
with our conclusions. 

The variance in the number density from the mocks is 
11%, somewhat better t han the 18% uncertainty on «drg from 
iMarchesini et al.l [20071 owing to the larger sample size. We 
use this 11% uncertainty when modeling the clustering. 

3. RESULTS 
3.1. Fitting the Quadri et. al. data 

We use a Monte Carlo Markov Chain (MCMC) analysis to 
determine the best-fit model parameters as well as their uncer- 
tainties. Using w{9) and «drg, the model with the minimum 
has parameters /^sat = 0.69, /smax = 0.99, and k = 0.48, 
yielding a number density of 6.4 x 10""* and /icen = 1-46. The 
for this model is 7.0 with 13-3 degrees of freedom. Fig- 
ure[3]compares our best-fit model with the Q08 data. At most 
scales the model is within the I- a errors. The highly cor- 
related nature of the data in the two-halo regime reduces the 
significance of the offset between the model and the data. The 
thin dotted lines show the two models presented in Q08. The 
upper curve, which matches the data at large scales, diverges 
rapidly from the data at < 100 arcsec because nearly all 
galaxies in this model are satellites in high-mass halos. The 
lower curve matches the data adequately on small scales, but 
is significantly below the data in the two-halo term, even with 
our new errors estimates. 

The differences between these two analyses are driven by 
the treatment of the second moment at small scales and the 
choice of halo bias models at lar ge scales. The calcul ations 
of Q08 are based on the model of iHamana et al.l (l2004l) (also 
used by lLee et al.ll2006l to model the cluster ing of LBGs, and 
similar to the model of iBullock et al.ll2002l) . In this formula- 
tion, the mean occupation function is a simple power law of 
the form 



_ J {M/Mif if M>M, 



if M< M„ 



(6) 



Even though equation (|6]l has one fewer free parameter than 
our {N)m in §2.2, it has significantly less freedom because 
it does not parameterize centrals and satellites separately. It 
is well known that assuming entirely Poisson fluctuations 
around {N)m produces clustering results that cannot match 
observations; the second moment is too high at low {N)m 



(see, e g..lBenson et al.ll2000l:ISeliakll2000l:IScoccimarro et al.1 
l200lt iBerlind & Weiiibergll2002h "n the HOD model of §2, 
the sub-Poisson scatter is a natural consequence of parameter- 
izing the centrals and satellites separately; the second moment 
of the satellites is always Poisson, but when (A^sat)M falls be- 
low unity the total distribution of pairs becomes sub-Poisson 
because there is no scatter ( or only a nearest-integer distribu- 
tion) of central galaxies (see lKravtsov et al.l2004tlZheng et all 
l2005h . In addition, the separation distribution of central- 
satellite pairs is different than satellite-satellite pairs. In the 
former, the pair distribution follows the halo density profile. 
In the latter, the pair distribution is represented by the den- 
sity profile convolved with itself. In galaxy distributions with 
high satellite fractions (say, ^ 30% for L < galaxies), the 
overall shape of the one-halo term will differ from samples 
with a low fraction of satellites (such as LRGs, DRGs, or any 
other commensurately bright sample that populates the high- 
mass end of th e halo mass functio n; see the discussion in the 
Appendix A of lZheng et ani2008h . 

In equation ^ there is no explicit separation of central and 
satellite galaxies. Thus, sub-Poisson fluctuations are speci- 
fied in the model ad- hoc when {N)m < 1- More importantly, 
the model specifies a binary switch between radial pair distri- 
butions; at {N)m > 1, the entire one-halo term is calculated 
assuming a satellite-satellite profile even though, at that mass 
scale, almost all pairs are still central-satellite. In Figure [3] 
the relative contribution of central-satellite pairs, S,i'l{r), and 
satellite-satellite pairs, Cih('")' to the one-halo term is shown 
with the dashed and dot-dashed curves, respectively. At all 
scales Cihi'') dominates because the overall satellite fraction 
of this sample is low (16%); although satellites dominate the 
pair counts at M > 2M[, the number of halos at that mass is 
nearly negligible because Mi is on the exponentially-falling 
tail of the mass function for this samp le. One can see that 
the shape of the one-halo term from the iHamana et al.l (l2004h 
model is the same as ^ih('"). This gives the one halo term 
a higher amplitude and larger radial extent than a model in 
which the two one-halo terms are properly weighted. Thus, to 
match at small scales, the best-fit model forces the large-scale 
clustering lower than it otherwise would be. 

The differences between Q08 and this work at small scales 
do not entirely account for the differences at large scales. Half 
the difference can be accounted for by t he choice o f large- 
scale bias. Q08 use the bias model of ISheth et all (l200Ih 
(hereafter SMT), which was calibrated on small-volume sim- 
ulations. There is debate in the literature over the accuracy of 
the SMT bias function at high masses (relative to the non- 
linear mass scale M^.) compared to the spherical collapse 
model (e.g.. lCole & Kaiseii[l989l: iMo & White|[T996l). which 
predicts stronger bias at these mass scales (see lCohn & White! 
[2008; Reed et al. 2008). Our model is intermediate between 
these two models and agrees perfectly w i th the bias of the Mil- 
lennium simulation shown in !Gao et all (!2005!) . Our new bias 
prescription increases the large-scale bias of the galaxy sam- 
ple by 10% relative to SMT, increasing the large-scale ampli- 
tude of ^(r) by 20%. The results of Figure [TJ) also show that 
our bias prescription is adequately describing the clustering of 
halos in the L120 simulation. Using the SMT bias function, 
our best-fit model has a = 10.0, while the Q08 model yields 
X^ = 26 A. Using the SMT bias function has no effect on the 
parameter constraints in §3.3. For comparison, using diagonal 
errors only results in x^ values of 10.9 for the best-fit model, 
19.5 when using SMT bias, and 32.0 for the Q08 model. 
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Fig. 4. — Two examples of the combined constraints of w{9} and 3>drg- Tlie top row sliows results from assuming tiiat DRGs are a random sample of all 
galaxies. The top left panel shows the HOD; the solid curve shows (Aten)** for DRGs, and the dashed line shows (A'sat)Af for DRGs. The dotted curve shows 
(N)m for all galaxies for comparison. The middle panel shows the luminosity function of Marchesini et al. (20()^ (points with errors). The solid hne is the model 
luminosity function. The dotted line is the luminosity function for all galaxies. The thick solid line is the completeness limit of the Q08 data. The right panel is 
the clustering of DRGs. The points represent the Q08 data, while the curve is the model prediction from the HOD in the far left panel. Bottom row: The best-fit 
model from the combined constraints of w(9) and "^drg ■ The mean mass of central DRGs is 2 1 % higher than the overall sample of central galaxies (~ 40% 
higher than non-DRG central galaxies). This model produces a better fit to both the luminosity function and w(8). 




Fig. 5. — Panel (a): Constraints on /rsm and /icen from the MCMC chains. Thick contours represent the I— a and 2 — rT constraints using both w^B) and 
'J'drg- Thin contours correspond to using w(6) alone. The gray dot represents the value /rsm obtained if equation j4) were applied to the subhalos. The white 
dot represents a model in which DRGs are a random sample of all galaxies. Panel (b): Constraints on /^sat and /smax. Contours are as in panel (a). The dotted 
line indicates the overall fraction of DRGs, 44%. Panel (c) shows the constraints on //jsat and (/ccn), the mean central DRG fraction integrating over all halos. 
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3.2. Including the Luminosity Function 

To a first approximation, rest-frame 7?-band absolute mag- 
nitude should be monotonically correlated with dark matter 
halo an d subhalo mass. Th is is borne out by the SHAM re- 
sults of lConrov et al.l (l2006l) . and it is further supported by the 
comparison betwe en SHAM and the c lustering of z ^ 2 star- 
forming galaxies jConrov et al.ll2008l) . By construction, the 
luminosity function of the galaxies in our mo del exactly fol- 
low that measured bv lMarchesini et aT1(l2007h for all galaxies 
(DRG and non-DRG). For each galactic halo, the luminos- 
ity of the galaxy within it is set by matching the cumulative 
halo abundance to the cumulative luminosity function. We 
can combine the luminosity set in this way with the mass- 
dependent DRG occupation function described in §2.2 to cre- 
ate a DRG luminosity funct ion, which has also been measured 
bv lMarchesini et aTl (|2007|) . To do this, we use the halos and 
subhalos identified in the LI 60 simulation. Each halo has a lu- 
minosity obtained through the abundance-matching method. 
Whether or not the galaxy is a DRG is decided by Monte 
Carlo based on /«cenW and /ssat- 

Figure |4] shows two examples of the combined constraints 
of w{0), «DRG and <1>drg^- In the top row, we show a model 
in which DRGs are simply a random sample of all galaxies. 
Thus /^cen = 1 (k = 0) and /^max = /j?sat = 0.44. The occupa- 
tion function for DRGs is the same as {N)m for all galaxies 
but shifted down by 0.44. The luminosity function produced 
by this model is in good agreement with the observations, but 
the clustering is clearly inconsistent with the measurements, 
both at large and small scales. The bottom row presents the 
results from our best-fit model from this combined approach. 
In this model, k = 0.18, /^sat = 0.91 and /^max = 0.66, yield- 
ing /icen = 1.21. This shift of the mass scale of central DRGs 
produces a noticeable increase in the large-scale bias of w{d). 
The increased fraction of red satellites increases the amplitude 
of one-halo term to agree well with the data. The luminosity 
function of this model agrees essentially perfectly with obser- 
vations, producing the right number of high-luminosity ob- 
jects while still matching the faint end down to the complete- 
ness limit of the Q08 sample. With 17-3 degrees of freedom, 
this model yields x" lv = 1 1 jW. 

3.3. Parameter Constraints 

We can use the observational data to place constraints on k, 
/ssat, and fumax- Figure |5] shows the results from the MCMC 
chains. The thick contours are the 1- and 2-cr constraints 
using 1^(6*), «DRG, and $drg, while the thin contours repre- 
sent results using only w{9) and «drg- Figure plots the 
constraints on the /fisat-Mcen plane. It is clear from Figure |4] 
that a large fraction of satellites must be red in order to match 
the clustering; the best-fit value of /^sat = 0.9 and the 2-cr 
lower limit is ~ 50% when marginalizing over other parame- 
ters. Ratios of the halo mass scale between central DRGs and 
all central galaxies above /^tcen =~ 1-5 are excluded at > 2-(t, 
but a model in which /icen = 1 (i e-, central DRGs are a random 
subsample of all central galaxies) is also excluded at roughly 
2-cr. 

FigurelSj) shows the parameter constraints in the /ssat- /«max 
plane. The dotted line indicates the mean DRG fraction. This 
figure highlights the benefit of using the luminosity function 
as an additional constraint on the models. When considering 

' We include both jiDRG and 'I'drg as independent data because <&drg 
does not necessarily include noRO- We find little difference in parameter 
constraints when excluding "drg- 
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Fig. 6. — Cumulative distribution of subhalos with accretion times (time 
since accretion) longer than t. Histogram shows the numerical data. The 
dotted line is a fitting function of the form exp(— 3.90- The shaded region 
shows the 2-a constraints on the fraction of satellites that are required to 
be DRGs to match the data. The thick solid line is the best-fit value; dashed 
lines represent the 1- and 2-a confidence levels after marginalized over other 
parameters. 



w(9) and «drg only, the model strives to push the large-scale 
bias as high as possible. This results in a best-fit /sniax of 
unity, implying all halos above M ~ lO'^'^ /i~'M0 contain a 
DRG at their center However, this overproduces the num- 
ber of bright DRGs at the expense of the low-luminosity end. 
Using w{9) and $drg, the best-fit value of /^max is 0.66. 

Figure|5}; shows the constraints on red satellite fraction and 
the mean red central fraction averaged over all halos. The 
constraints in all three panels imply that a higher fraction of 
satellites than centrals must be red in order to match the data. 
(We remind the reader that /«max is the asymptotic value of 
fuceniMy, at fixed /ssat, the higher the value of fumax, the lower 
the red fraction at lower halo masses in order to keep the num- 
ber of central DRGs constant. The overall red central fraction 
will usually be near 44% because the contribution of satel- 
lites is small.) The gray circle in Figures|5^ and|5]3 shows the 
/ssat that results if we apply the mass-dependent /scenC-^) to 
the subhalos themselves for the best-fit value of k. Because 
subhalos have a lower mean mass than parent halos, in this 
case the fraction of red satellites is less than the fraction of 
red centrals if they have the same duty cycle. Thus, using 
this model, the fraction of red satellites is 32%, a value that is 
clearly outside the 2-cr contour. 

4. DISCUSSION 

4.1. Satellite DRGs 

The results from Figure |5] imply that the physical processes 
that turn a galaxy red, either by creation of dust or through star 
formation history, are different for centrals and satellites at 
z > 2. At low redshift such a scenario is expected; after accre- 
tion, gas is ram-pressure stripped and star formation in satel- 
lites is "quenched". However, the situation at high redshift is 
quite different; the universe at z = 2.3 is only ~ 2.9 Gyr old in 
our cosmology, and most mergers of small halos onto bigger 
ones have occurred quite recently. Thus, one would naively 
expect that red satellites were red before they were accreted 
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the Q0 8 sample. Large circles and their eiTors aie taken from lOuadri et all 
120071) . which show the clustering of J-K < 2.3 galaxies from a smaller 
sample. The lines are HOD predictions for the clustering of non-DRGs from 
the best-fit value of fi^m but with varying red satellite fractions. The dotted 
line is the prediction of the overall best-fit model (bottom row of Figure|4]) 
The small squares are the measurement from Ouadri et al. (2007) for DRGs, 
which are consistent with the Q08 data (black circles connected by the dash- 
dot curve). 



(or would have become red regardless of accretion); i.e., that 
the processes that make a galaxy red are driven primarily by 
the mass of its galactic halo rather than their position as an 
accreted object. Such a scenario is difficult to reconcile with 
the Q08 data. 

Figure |6] shows the cumulative distribution of satellite ages 
(i.e., time since they were accreted) for subhalos above the 
minimum mass threshold of (A^) ^ in the simulations. Because 
the accretion times are discretized, the dotted line is a fitting 
function of the form exp(-3.9f) that allows for quick interpo- 
lation between the points. Naively we may infer a quenching 
time scale based on the fraction of satellites being DRGs (e.g., 
200 Myr from the l-a line in Fig. 6) , but this makes the strong 
assumption that all satellite DRGs have little to no star forma- 
tion. Let us assume that 3 0% of D RGs have significantly at- 
tenuated star formation tLabbe et al..2005 '), and 32% of DRG 
satellites were classifiably DRGs before accretion (from the 
gray dot in Figure |5j. Using the best-fit value of /^at = 0.9, 
0.9 X (1 -0.32) X 0.3 = 0.18 is the fraction of all galaxies that 
have their star formation attenuated through satellite accre- 
tion. From Figure |6] a star formation quenching timescales 
of 450 Myr is obtained. This is is clear contrast with re- 
sults from low redshift, where the exponential timescale for 
the r eduction of star formation in satellite s is closer to 2-2.5 
Gyr jWang et al.ll2006t iKimm et al.ll2008l) . These results can 
be brought into better agreement if the exponential timescale 
scales with the dynamical time of the galactic halo, reducing 
the e-folding time by (1 -l-z)'^ = 6, depending on the number 
of e-folds required to be considered "red and dead". Robust 
conclusions are difficult to establish given the uncertainty in 
the fraction of DRGs with no star formation, uncertainty in 
how that fraction is distributed between centrals and satellites, 
and the uncertainty in /ssat- 

Regardless of current star formation rates, it is also difficult 



to reconcile the high overall fraction of satellite DRGs with 
satellite accretion timescales. While there is a physical mech- 
anism that halts star formation for satellite galaxies, accretion 
onto a larger halo cannot account for excess A'-band emission 
from an advanced (multi-Gyr) stellar population that is part of 
the definition of a DRG. Irrespective of accretion time, the old 
stars were formed before the galaxy became a satellite galaxy. 
Thus the satellite had to know it would be accreted well be- 
forehand. For z = 2.3 10'^ h~^Mp) halos, th eir average mass 
was 1/10 at z 5 dWechsler et al.ll2002l) . but their rapid 
growth and high relative mass can influence their larger envi- 
ronment and t he smaller halos that will eventually be accreted 
at later times (I Wang et al.ll2007l: iHahn et alJl2008HD"alal et alJ 
I2008h . These lower- mass halos within the lagrangian radius 
of the large halo will have formed earlier and possibly have 
older stellar stellar populations than a halo of the same mass 
which formed outside the lagrangian radius of a large halo. 
By z w 2.3, these early-forming halos have become satellites 
in the nearby high-mass halos. 

It would be perhaps too coincidental that this effect would 
only alter the formation of halos that would become subhalos 
at the epoch of the Q08 sample, especially given that most ac- 
cretion events are very recent. High mass halos at z = 2.3 will 
continue to rapidly grow and alter the formation trajectories of 
smaller halos around them, inducing assembly bias — the cor- 
relation between halo occupation and large-scale environment 
at fixed halo mass — cited by Q08 to explain the high cluster- 
ing amplitude of DRGs. If this assembly bias extends to halos 
outside the virial radii of M > 10^^ Mq halos (but within 
their vicinity), the large-scale clustering of DRGs would be 
enhanced, supporting the conclusions of Q08. The data do 
not require any assembly bias to fit the large scale w{9), but 
it is implied by the high fraction of satellite galaxies being of 
DRGs. 

Recent developments in stellar population synthesis, specif- 
ically the inclusion of thermally-pulsating AGB stars, may al- 
leviate some of the tension between the age of the universe at 
z = 2.3 and the stellar populati on ages required to reproduce 
the observed colors of DRGs dTonini et alj|2008l) . although 
this mechanism is not related to the high fraction of satellites 
being DRGs or the fraction of DRGs without star formation. 

4.2. Central Galaxies 

For central galaxies, the luminosity function alone supports 
a model in which DRGs are a random subsample of the ha- 
los occupied by bright galaxies. In contrast, the clustering of 
DRGs supports a model in which massive halos (M > 10'^'^) 
only house DRGs at their centers. Considered together, a 
more intermediate picture becomes clear in which massive ha- 
los host DRGs ~ 2/3 of the time, in good agreement with the 
mass-selected sample of 'van Dokkum et al. (2006). We have 
made two strong assumptions when predicting $drg from our 
models: 1) that there is no scatter in the relation between halo 
mass and luminosity, and 2) that mass and 7? -band magnitude 
are monotonically correlated regardless of galaxy color. For 
1), we note that scatter does not change the predicted $drg 
when /icen and /^max are high; these models are excluded when 
using $DRG (cf, Figure|5j5)- In such models all high-mass ha- 
los contain DRGs at their centers, thus scatter only changes 
what halo each DRG lives in and not the distribution of DRG 
luminosities on the bright end. Significant scatter between 
halo mass and luminosity at low luminosities would decrease 
the observed amplitude of the correlation function, contrary to 
the observations (unless assembly bias is canceling out this ef- 
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feet). For 2), with the given data we cannot rale out a model in 
which red galaxies and blue galaxies of the same R magnitude 
occupy halos with substantially different masses. Our model 
makes a clear prediction for the clustering of blue galaxies 
that is sensitive to this assumption. For a complete model, the 
clustering of blue and red galaxies should be modeled simul- 
taneously. 

Fi gure [7] shows the clu stering of non-DRGs (J-K < 2.3) 
from lOuadri et al.l (l2007h . a distinct and much smaller sample 
of galaxies than Q08. The three lines are predictions from our 
best-fit value of k = 0. 18 for different values of /ssat- At large 
scales, the measured amplitude of w{6) is consistent with our 
model predictions. This agreement supports the second as- 
sumption above and argues against significant assembly bias 
in the DRG population; if the clustering of DRGs is enhanced 
by assembly bias then the clustering on non-DRGs must be 
suppressed by the same effect. At small scales, there is a sig- 
nature of one-halo clustering in the blue galaxies. A model 
with nearly no blue satellites appears to be difficult to recon- 
cile with the observations, but these data cannot distinguish 
between models in which 90% and 50% of the satellites are 
red. Because this is a different galaxy sample from Q08 with 
a different redshift distribution, the comparison in Figure |7] 
is meant to be qualit ative only. We note , however, that the 
DRG clustering from lOuadri et al.l (l2007h is consistent with 
that measured in Q08 (dash-dot curve and squares in Figure 
|7]). Measurements of the clustering of blue galaxies within the 
Q08 sample will enhanced our constraints both for halo occu- 
pation and for assembly bias. The relative bias between red 
and blue galaxies within the same field can also shed light on 
the issue of assembly bias. Because they would be measured 
from the same volume, large-scale modes would affect the 
clustering of both red and blue galaxies in the same way and 
the relative b i as is therefore independen t of cosmic variance 
(ISeliakll2008l:lMcDonald & Seliakll2008h . Assembly bias, on 
the other hand, would increase the clustering difference be- 
tween red and blue subsamples. 

5. CONCLUSIONS 

The clustering of DRGs measured by Q08 can be ade- 
quately fit by standard halo occupation models. With proper 
estimates of the cosmic variance at large scales and a more 
robust implementation of the HOD, our best-fitting model 
has x"/^ < 1- Thus, from the large-scale clustering alone 
there is no compelling evidence for assembly bias in the 
halo occupation of DRGs. The high clust ering amplitude of 
DRGs has been mea sured by other author s dDaddi et al.ll2003l : 
iGrazian et ani2006h . but the iDaddi etlTI (l2003h data is fully 
described by HOD modeling as well dZhengl2004l) . Given that 
the Q08 sample is by far the largest at present, the inclusion 
of other currently available data sets is not likely to change 
our conclusions, though cosmic variance is still an important 
systematic on large scales. 



To fit the prominent break in the w(9) between the two- 
halo and one-halo terms, a large fraction of satellite galax- 
ies above the relevant luminosity threshold are required to be 
DRGs. The best-fit model has a red satellite fraction, /ssat, of 
nearly unity, with a 2-0- lower limit of 50%. This implies that 
the mechanisms through which DRGs form are more efficient 
(or more frequent) for satellite DRGs than for centrals, even 
though the accretion times of most subhalos are very recent at 
these epochs. At the best-fit value of /«sat = 0.9, the timescale 
for star-formation quenching is roughly 450 Myr if quenching 
begins at the accretion time. Regardless of the star formation 
rates at z = 2.3, to produce the evolved stellar populations in 
satellite DRGs, subhalos are required to know that they will 
become subhalos well before they are accreted. This implies 
these object have some knowledge of the large-scale environ- 
ment beyond that of their host halo mass (i.e. assembly bias), 
but further data is required to quantify the effect. 

Down to the completeness limit of the Q08 data, central 
DRGs have an average halo mass that is 25% higher than the 
overall sample of galaxies (thus 50% higher than non-DRGs), 
but at the highest halo masses DRGs do not occupy all halos; 
~ 1/3 of the most massive halos at z = 2.3 still contain blue 
galaxies at their centers. 

Clustering data for high-redshift galaxies are becoming 
precise enough that realistic models for the occupation of 
galaxies within halos are important for robust interpretation. 
The approach taken in this work for modeling halo occu- 
pation is more realistic than earlier mod els that were devel- 
oped to modeling hi gh-redshift data (e.g. iBullock et al.ll2002l : 
iHamana et al.ll2004l) . Explicit treatment of central and satel- 
lite galaxies is essential in modeling the transition regime be- 
tween one and two halo pairs, which is particularly important 
as a constraining feature now that it is well measured in the 
data. The approach taken here, which combines information 
from simulations which resolve substructure to constrain the 
global sample with a flexible approach for modeling color de- 
pendence, is particularly powerful, and can also be applied to 
other color-selected galaxy samples. 
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APPENDIX 
A. HALO BIAS 

The large-scale bias of halos used in this paper is based on the spherical overdensity halo catalogs of [Tinker et al.l ( l2008al) . 
The full results fo r both large-scale bias and scale-dependent bias as a function of halo overdensity will be given elsewhere 
dTinker et alj|2009l) . In the Tinker et. al. results, the fitting function for large-scale bias is 

b(i^)=l-A^^+Biy'' + Ci^'. (Al) 
where ly = Sc/(7(M), 6e = 1.686 and <t is the linear matter fluctuations within the lagrangian radius of a halo of mass M. For 
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A = 200 halos, A = 1 .04, a = 0. 132, = 0. 183, ^? = 1 .5, C = 0.262, and c = 2.4. 
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